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ABSTRACT 

Next generation architectures necessitate a shift away from 
traditional workflows in which the simulation state is saved 
at prescribed frequencies for post-processing analysis. While 
the need to shift to in situ workflows has been acknowledged 
for some time, much of the current research is focused on 
static workflows , where the analysis that would have been 
done as a post-process is performed concurrently with the 
simulation at user-prescribed frequencies. Recently, research 
efforts are striving to enable adaptive workflows, in which the 
frequency, composition, and execution of computational and 
data manipulation steps dynamically depend on the state of 
the simulation. Adapting the workflow to the state of sim¬ 
ulation in such a data-driven fashion puts extremely strict 
efficiency requirements on the analysis capabilities that are 
used to identify the transitions in the workflow. In this pa¬ 
per we build upon earlier work on trigger detection using 
sublinear techniques to drive adaptive workflows. Here we 
propose a methodology to detect the time when sudden heat 
release occurs in simulations of turbulent combustion. Our 
proposed method provides an alternative metric that can be 
used along with our former metric to increase the robustness 
of trigger detection. We show the effectiveness of our metric 
empirically for predicting heat release for two use cases. 

1. INTRODUCTION 

Concurrent analysis frameworks have been developed to 
process raw simulation output as it is computed, decoupling 
the analysis from I/O. Operations sharing primary resources 
of the simulation are considered in situ, while in transit pro¬ 
cessing involves asynchronous data transfers to secondary 
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resources. Both in situ |15| [6] [8] and in transit [14} jTJ [2] 
workflows perform analyses as the simulation is run and pro¬ 
duce results, which are typically much smaller than the raw 
data, mitigating the effects of limited disk bandwidth and 
capacity. Concurrent analyses are often performed at pre¬ 
specified frequencies, which is viable for analyses that are 
not too expensive - in terms of runtime (with respect to a 
simulation time step), memory footprint, and output size. 
However, for many other analyses that are resource intense, 
prescribed frequencies will not suffice because the scientific 
phenomenon being simulated typically does not behave lin¬ 
early (e.g., combustion, climate, astrophysics). When the 
prescribed I/O or analysis frequency is high enough to cap¬ 
ture the events of interest, the costs incurred are too great. 
On the other hand, a cost-effective, lower analysis-frequency 
may miss the scientific events that simulation is intended to 
capture. 

An alternative approach is to perform expensive analy¬ 
ses (denoted A e ) and I/O in an adaptive fashion, driven 
by the data itself. A domain-agnostic approach is presented 
in |l2 Hi, based on entropy of information in the data. Our 
previous work in 3 provides a framework for making data- 
driven control-flow decisions that can leverage the scientists’ 
intuitions, even when the algorithms to capture those intu¬ 
itions would otherwise be too expensive to compute. Our 
methodology involves a user-defined indicator function that 
is computed and measured in situ at a relatively regular and 
high-frequency (i.e., greater than the frequency with which 
the I/O or A e would be prescribed). Along with the indi¬ 
cator, the application scientist defines an associated trigger, 
a function that returns a boolean value indicating that the 
indicator has met some property, for example a threshold 
value. While our methodology is intuitive and conceptually 
quite simple, the challenges lie in defining indicators and 
triggers that capture the appropriate scientific information, 
while remaining cost efficient in terms of runtime, memory 
footprint, and I/O requirements, so that they can be de¬ 
ployed at the high frequency that is required. 

In 3 we show that chemical explosive mode analysis, de¬ 
noted CEMA, can be used to devise a noise-tolerant in¬ 
dicator for heat release (the quantity of interest that the 
scientists would like to capture), thus making it a good 
candidate to drive adaptive workflows. However, exhaus¬ 
tive computation of CEMA values dominates the simula¬ 
tion time. To overcome this bottleneck, we proposed a 
quantile sampling approach with provable error/confidence 
bounds. These bounds depend only on the number of sam¬ 
ples and are independent of the problem size. We also de- 



signed an indicator, referred to as the P-indicator, based on 
the quantile sampling approach. Our experiments on homo¬ 
geneous charge compression ignition (HCCI) and reactivity 
controlled compression ignition (RCCI) simulations showed 
that the proposed method can detect rapid increases in heat 
release and is computational efficient. 

In this paper, we propose an alternative indicator, referred 
to as the C-indicator, and associated trigger function. The 
proposed technique is based on the coefficient of variation 
and similar in essence to our earlier technique that used 
quantiles, as it tries to detect shrinking in the range of up¬ 
per percentiles in the distribution of CEMA values. The 
coefficient of variation among the top quantiles decreases as 
the range of the values shrink and can serve as a detector 
for the subsequent heat release that the scientists wish to 
capture. The new technique based on coefficient of varia¬ 
tion is not proposed as an alternative to replace our former 
metric. Rather, we believe that a set of trigger detection 
mechanisms, when used collectively, can provide more ro¬ 
bust detections, as collectively they will be less prone to 
false positives. Our experiments on HCCI and RCCI simu¬ 
lations show that the proposed method can efficiently detect 
rapid increases in heat release. 

2. COMBUSTION AS A USE CASE FOR 
TRIGGER DETECTION 

We demonstrate our approach applied to a combustion use 
case, using S3D f], a direct numerical simulation of combus¬ 
tion in turbulence. The combustion simulations in our use 
case pertain to a class of internal combustion (IC) engine 
concept. We are interested in two specific techniques called 
homogeneous-charge compression ignition (HCCI) 4] and 
reactivity-controlled compression ignition (RCCI) [9 5 . In 
both cases, heat release starts in the form of small kernels 
at arbitrary locations. Eventually, multiple kernels ignite 
as the overall heat release reaches a global maximum and 
subsequently declines. Since these simulations are computa¬ 
tion and storage-intense, we want to run the simulation at 
a coarser grid resolution and save data less frequently dur¬ 
ing the early build-up phase. When the heat release events 
occur, we want to run the simulation at the finest grid gran¬ 
ularity possible, and store the data as frequently as possible. 
Therefore, it is imperative to be able to predict the start of 
the heat release event using an indicator and trigger that 
serve to inform the application to adjust its grid resolution 
and I/O frequency accordingly, see Figure [l] 

3. DESIGNING A NOISE-RESISTANT 
INDICATOR AND TRIGGER 

In this section we describe the design of an indicator and 
trigger for heat release for our combustion use case. To 
provide context, we begin by discussing the intuitions that 
informed our design. 

3.1 Chemical Explosive Mode Analysis 

Chemical explosive mode analysis (CEMA) is known to be 
a reliable technique to predict incipient heat release. Here 
we provide brief description and refer to 10, 13 for de¬ 
tails. The conservation equations for reacting species can 
be written as ^ = a>(y) + s(y). The vector y represents 
temperature and reacting species mass fractions, uj is the 
reaction source term and s is the mixing term. Then, the 



Figure 1: The minimum (blue), maximum (red) and mean 
(black) heat release values for each time step in the simulation. 
Early in the simulation, we want to run at a coarser grid and 
save data less frequently. When heat release occurs, we want to 
the finest grid and save the data as frequently as possible. The 
vertical lines define a range of time steps within which we would 
like to make this transition, as identified by a domain expert. 

Jacobian is + J s where = d ^' > and J s = ds ^ . 
The eigen-decomposition of the chemical Jacobian, can 
be used to infer chemical properties of the mixture. Let A e 
be the eigenvalue with the largest real part. A e is defined 
as a chemical explosive mode (CEM) if Re(A e ) > 0, which 
indicates that point will undergo ignition. If it has under¬ 
gone ignition, we have Re( A e ) < 0. The presence of a CEM 
indicates the propensity of a mixture to ignite. 

Our CEMA-based indicator is based on global trends of 
CEMA over time. Consider Figure [2] which provides a sum¬ 
mary of the trends of CEMA values across all time steps in 
a simulation. At timestep t, let C(t) be the sorted (in non¬ 
decreasing order) array of CEMA values on the underlying 
mesh. For a £ (0,1], the a-percentile is the entry C(t) |- a jvi ■ 
More specifically, it is the value in C(t) that is greater than 
at least [aiV] values in C(t). We denote this value by p a (t). 



Figure 2: Percentile plot of CEMA values. The blue curves 
correspond to po.Ol,0,02,...,o.lj the green curves to po.2,o.3,...,0.9 
and the red curves to po. 91 , 0 . 92 ,....l- 

We notice that as the simulation progresses, the distance 
between the higher percentiles (red curves) decreases then 
suddenly increases. This is illustrated in the plot by the 
spread in the red curves that occurs between the dashed 
lines defined by the domain expert that indicate the desired 
transition phase. Our goal is to capture this empirical ob¬ 
servation with a metric, and use it to predict heat release. 

























Our empirical observation is consistent with the underly¬ 
ing physics, and we refer to [3] for a detailed discussion on 
the underlying physics. 

3.2 Designing a Noise-Resistant CEMA-Based 
Indicator 

We introduce an indicator function, C-indicator that quan¬ 
tifies the distribution of the top quantiles of CEMA values 
over time. The C-indicator is equal to the coefficient of vari¬ 
ation (COV) of the top percentiles of CEMA values, i.e. the 
ratio of their standard deviation to their mean. 

Formally, quantiles are defined as values taken at regu¬ 
lar intervals from the inverse of the cumulative distribution 
function of a random variable. For a given data set, quantiles 
are used to divide the data into equal sized sets after sorting, 
and the quantiles are the values on the boundary between 
consecutive subsets. A special case is dividing into 100 equal 
groups, when we can refer to quantiles as percentiles. This 
paper focuses on percentiles with numbers in the [0,1] range 
(although all techniques presented here can be generalized 
for any quantiles). For example, the 0.5 percentile will refer 
to the median of the data set. 

We define our indicator using 2 parameters: a and /3. 
We want to detect whether the range covered by top quan¬ 
tiles shrinks, and a represents the lower end of the top 
percentiles, whereas /3 represents the top percentile consid¬ 
ered. Therefore, the range of top percentiles we measure is 
\p a (t),pp(t)\. In our indicator, we choose a < /3 (typically 
in the range [0.90,0.99]). We measure the spread at time t 
by the COV of the p s {t) values for a < s < /3 which results 
in the following C-indicator: 
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where 

( 2 ) 

S = Oi 

and N is the number of percentile values considered between 
a and /3. For instance, if a = 0.90 and /3 = 0.99, we would 
obtain N = 10. 

Figure [3] illustrates percentile plots for heat release (top 
row) and CEMA (middle row). In the percentile plots, the 
lowest blue curve and the highest red curve correspond to 
the 1 and 100 percentiles (po.oi and pi), respectively. The 
blue curves correspond to po.oi,o. 02 ,...,o.i, the green curves 
to Po. 2 ,o. 3 ,...,o .9 and the red curves to po. 9 i,o. 92 ,...,i- The C- 
indicator evaluated for a = 0.92 and /3 = 0.99 is shown in 
the bottom row of Figure [3] Results are generated for four 
test cases described in Table [T] The vertical dotted lines 
were identified by a domain expert who, via examination 
of heat release and CEMA percentile plots, visually located 
the range of time steps in the simulation where the mesh 
resolution and I/O frequency should be increased. We refer 
to this range of time steps as the “true” trigger range we 
wish to identify with the C-indicator and trigger functions. 
Note, for the RCCI cases, there are two ignition ranges. To 
simplify the following exposition, we focus on the second 
rise in the heat release rate profiles, as this is the ignition 
stage of interest to the scientists. However, we note that our 
approach is robust in identifying the first ignition stage as 
well. 
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Figure 3 : Plots showing the (top row) percentiles of the heat 
release rate, (middle row) percentiles of CEMA, and (bottom 
row) the C-indicator, as indicated. The blue curves correspond 
to po.01,0.02,...,0.1, the green curves to po.2,0.3,...,0.9 and the red 
curves to po.9i,o.92,...,i• The C-indicator shown is evaluated for 
ao = 0.92 and /3 = 0 . 99 . The vertical dotted lines crossing the 
images indicate a window of acceptable “true” trigger time steps, 
as identified by a domain expert. For the RCCI cases, the trig¬ 
ger time ranges are based on the High Temperature Heat Release 
(HTHR), i.e., the second peak in the Heat Release Rate (HRR) 
profiles. 


3.3 Defining a Trigger 

In addition to defining a noise-resistant indicator func¬ 
tion, we also need to define a trigger function that returns a 
boolean value, capturing whether a property of the indicator 
has been met. Looking at Figure [3] we notice that across 
all experiments from Table [I] the C-indicator is increasing 
during the true trigger time step windows. Therefore, we 
seek to find a value re £ (0,1), such that C a ,p{t) crosses rc 
from below, as the simulation time t progresses. 

Figure[l]plots the trigger time steps as a function of rc for 
Cc—o. 92 ,/ 3 —o. 99 (t). This plot shows that the viable range of 
rc is [0.01, 0.05] for the four use cases described in Table[l] 
The horizontal dashed lines indicate the true trigger range 
identified by our domain expert. We consider those values of 
rc that fall within the horizontal dashed lines to be good rc 
values for our trigger, with those values of rc falling below 
the dashed lines still considered viable. We note the plots for 
other values of a and /3 look similar and have been omitted 
due to lack of space in this text. 

4. COMPUTING INDICATORS AND 
TRIGGERS EFFICIENTLY: 

A SUBLINEAR APPROACH 

The previous section showed that C-indicator and trig¬ 
ger are robust to noise fluctuations and act as a precursor 
to rapid heat release in combustion simulations. In 3] we 
show the computational cost of computing CEMA-based in¬ 
dicators can be prohibitive for large-scale simulations (up 
to 60 times the cost of a simulation time-step depending 
on the simulation parameters). To mitigate the cost, we 
introduced a quantile sampling approach that comes with 
provable bounds on accuracy as a function of the number 
of samples. Most importantly, the required number of sam- 






















































Table 1: Four Combustion Use Cases analyzed in this study. The “true” trigger time ranges are estimated based on 95 — 100 th percentiles 
of the heat release rate. The computed time ranges were evaluated using our quantile sampling approach. The C-trigger and P-trigger 
were computed over 50 realizations of experiments with 20 samples per process. 


Problem 

Instance 

Number of 
Grid Points 

Number of 
Species 

“True” Trigger 
Time Range 

C-trigger 

Detection 

P-trigger 
Detection |5j 

Total 

Processes 

HCCI, T=15 

451,584 

28 

250-315 

235 - 265 

250 - 262 

1600 

HCCI, T=40 

451,584 

28 

175-225 

180 - 220 

213 - 220 

1600 

RCCI, case 6 

2,560K 

116 

38-50 

34 - 38 

28 - 45 

6400 

RCCI, case 7 

2,560K-10,240K 

116 

42-58 

32.5 - 35.5 

35 - 50 

6400 


a = 0.92, (3 = 0.99 



Figure 4 : This figure plots the trigger time steps as a function 
of T(j for C a= o.92,/3=o.99(0- There is a range of viable values of 
rc S [0.0, 0.05] that predict early stage heat release. 


pies for a specified accuracy is independent of the size of the 
problem, hence our sampling based algorithms offer excel¬ 
lent scalability. 

We summarize the quantile sampling method here and re¬ 
fer the reader to [3 for analytical and empirical proofs of 
convergence. Consider an array A € K^, in sorted order. 
Our aim is to estimate the o-percentile of A. (We use p a to 
denote the percentiles.) Note that this is exactly the entry 
A^ aN -\ ■ Here is a simple sampling procedure. 

1. Sample k independent, uniform indices r\, r 2 , ■ ■ ■, r*, in 
{1,2,..., N}. 

Denote by A the sorted array ]A(ri), A(r 2 ),..., A(rk)]- 

2. Output the a-percentile of A as the estimate, p a . 

We performed a series of experiments examining the vari¬ 
ation in the trigger time steps as a function of the number of 
samples used per process. The data for Figure [5] was gener¬ 
ated via 50 realizations of the C-indicator with a = 0.92 and 
/3 = 0.99, with rc drawn from [0.01,0.05]. The horizontal 
dashed lines in this figure define the range of true trigger 
time steps within which we would like to make the work- 
flow transition (as identified by a domain expert). This plot 
demonstrates that, even across the range of rc values, with 
a small number of samples per process, our quantile sam¬ 
pling approach can accurately estimate the true trigger time 
steps as defined by the domain expert. The accuracy of our 
quantile-based sampling predictions is also shown in Table[l] 
which lists the the triggers predicted by the C-indicator and 
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Figure 5: Plots illustrating the variability of the trigger time 
steps predicted by the C-indicator and trigger as a function of 
the number of samples per process. The data for these plots was 
generated via 50 realizations of the C-indicator with a = 0.92 and 
/3 = 0.99, and rc drawn from [0.01, 0.05].The horizontal dashed 
lines define the range of time steps within which we would like to 
make the workflow transition (as identified by a domain expert. 


P-indicator that we used in our previous work 3 . The P- 
indicator is computed as P a ,/3,-y = p^-p 1 2 , f° r ot — 0.94, 
p = 0.98, 7 = 0.01, and the trigger threshold chosen as 
[0.725,0.885]. The results show both methods are accurate. 
We are currently investigating how we can improve our re¬ 
sults by using two metrics concurrently. Note that using 
both metrics will not induce an additional burden, since the 
two metrics can be computed from the same set of samples. 

5. CONCLUSION 

We propose a new indicator and trigger for making data- 
driven control-flow decisions in situ. Using a provably robust 
sampling-based approach, of CEMA values (which, when 
computed in full, can cost up to 60 times a simulation time 
step). Our experiments show that our proposed indicator, 
based on the coefficient of variation, can efficiently predict 
rapid increases in heat release. We believe the new metric 
can be deployed collectively with a previously defined mech¬ 
anism, to provide more robust detections, since together 
they will be less prone to false positives. We note that using 
both metrics will not induce an additional burden, since the 
two metrics can be computed from the same set of samples. 
Our future work aims to explore deployment of both metrics 
jointly in production simulation runs. 
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